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Abstract 

We consider the problem of identifying patterns in a data set that 
exhibit anomalous behavior, often referred to as anomaly detection. In 
most anomaly detection algorithms, the dissimilarity between data sam- 
ples is calculated by a single criterion. When dissimilarities are calculated 
by multiple criteria, one might perform anomaly detection using a lin- 
ear combination of the multiple dissimilarities. If the importance of the 
different criteria are not known in advance, the algorithm may need to 
be executed multiple times with different choices of weights in the linear 
combination, perhaps selected by grid search. 

In this paper, we introduce a novel non-parametric multi-criteria ano- 
maly detection method using Pareto depth analysis (PDA). PDA uses the 
concept of Pareto optimality to detect anomalies under multiple criteria 
without having to run an algorithm multiple times with different choices 
of weights. The proposed PDA approach scales linearly in the number 
of criteria, unlike grid search methods, which scale exponentially. We 
find that PDA is able to outperform any linear combination of criteria, 
including weights selected by grid search. 



1 Introduction 

Anomaly detection is an important problem that has been studied in a variety 
of areas and used in diverse applications inc luding intrusion detection, fraud de- 



tection, and image processing (H odge and Austinl . 12004; .Chandola et al. . 200 



Many methods for anomaly detection have been developed using both para- 
metric and non-parametric approaches. Non-parametric approaches typically 
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involve the calculation of dissimilarities between data samples. For complex 
high-dimensional data, several dissimilarity measures corresponding to different 
criteria may be required to detect different types of anomalies. For example, 
consider the problem of detecting anomalous object trajectories in video se- 
quences. Multiple criteria, such as dissimilarity in object speeds or trajectory 
shapes, can be used to detect a greater range of anomalies than a single criterion. 
In order to perform anomaly detection using these multiple criteria, one could 
first combine the dissimilarities using a linear combination. However, in many 
applications, the importance of the different criteria are not known in advance. 
It is difficult to determine how much weight to assign to each dissimilarity mea- 
sure, so one may have to choose multiple weightings using, for example, a grid 
search. Furthermore, when the weights are changed, the anomaly detection 
algorithm needs to be re-executed using the new weights. 

In this paper we propose a novel non-parametric multi- criteria anomaly 
detection approach using Pareto depth analysis (PDA). PDA uses the concept 
of Pareto optimality to detect anomalies without having to choose weights for 
different criteria. Pareto optimality is the typical method for defining optimality 
when there may be multiple confiicting criteria for comparing items. An item 
is said to be Pareto-optimal if there does not exist another item that is better 
or equal in all of the criteria. An item that is Pareto-optimal is optimal in 
the usual sense under some combination, not necessarily linear, of the criteria. 
Hence, PDA is able to detect anomalies under multiple combinations of the 
criteria without explicitly forming these combinations. 

The PDA approach involves creating dyads corresponding to dissimilarities 
between pairs of data samples under all of the dissimilarity measures. Sets of 
Pareto-optimal dyads, called Pareto fronts, are then computed. Nominal and 
anomalous samples are located near different depths of Pareto fronts; thus us- 
ing the front depths of the dyads corresponding a test sample is an effective 
way to discriminate between nominal and anomalous samples. The proposed 
PDA approach scales linearly in the number of criteria, which is a significant im- 
provement compared to selecting multiple weights via a grid search, which scales 
exponentially in the number of criteria. Furthermore, we find that PDA outper- 
forms any linear combination of criteria when used with several state-of-the-art 
anomaly detection algorithms in two experiments involving both synthetic and 
real data sets. 

The rest of this paper is organized as follows. We first provide an introduc- 
tion to Pareto fronts and discuss related work in Section [2] We then introduce 
the concept of a dyad and relate Pareto fronts of dyads to the anomaly detection 
problem in Section [3l The PDA anomaly detection algorithm is then described 
in Section |3) Finally we present two experiments in Section [5] to evaluate the 
performance of PDA. 
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2 Background 
2.1 Pareto fronts 

The PDA method proposed in this paper utihzes the notion of Pareto optimal- 
ity, which has been studied in many apphcation areas in economics, computer 
science, and the social sciences among others. We provide an introduction to 
Pareto optimahty and define the notion of a Pareto front. 

Consider the fohowing problem: given n items, denoted by the set S, and K 
criteria for evaluating each item, denoted by functions fi, . . . , fx, select x & S 
that minimizes [/i(a:), . . . , /if (a;)]- In most settings, it is not possible to identify 
a single item x that simultaneously minimizes fi{x) for all i £ {!,..., K}. In 
order to decrease criterion i, i.e. by finding x* such that fi{x*) < fi{x), it may 
be necessary to increase criterion j, i.e. fj{x*) > fj{x). A minimizer can be 
found by combining the K criteria using a linear combination of the /^'s and 
finding the minimum of the combination. Different choices of weights in the 
linear combination could result in different minimizers; a set of items that are 
minimizers under some linear combination can then be created by using a grid 
search over the weights, for example. 

A more powerful approach involves finding the set of Pareto-optimal items. 
An item x is said to strictly dominate another item x* if x is no greater than x* 
in each criterion and x is less than x* in at least one criterion. This relation can 
be written as x >- x* if fi{x) < fi{x*) for each i e {!,..., K} and fi{x) < fi{x*) 
for some i. The Pareto front, which corresponds to the set of Pareto-optimal 
items, is the set of items in S that are not strictly dominated by another item in 
iS. It contains all of the minimizers that are found using linear combinations, but 
also includes other items that cannot be found by linear combinations. Denote 
the Pareto front by which we call the first Pareto front. The second Pareto 
front can be constructed by finding items that are not strictly dominated by any 
of the remaining items, which are members of the set S\J-i. More generally, 
define the ith Pareto front by 



For convenience, we say that a Pareto front Ti is deeper than J-j \i i > j. 

We illustrate the concept of Pareto fronts on a toy example. We draw 50 
samples from a uniform distribution on [0, 1] x [0, 1]. The ith sample is denoted 
by [xi,yi\. K — 2 criteria are used to evaluate each sample, where fi{[x, y]) = x 
and f2{[x,y]) — y. The first five Pareto fronts for one realization are shown 
in Figure [T] Notice that the fronts create layers of the data. Like an onion, 
the second layer is revealed by removing the first layer, and the third layer is 
revealed by removing the second layer, and so on. As we show in Section [31 the 
depths of the Pareto fronts carry valuable information for anomaly detection. 
Notice also that linear combinations can only identify four of the seven samples 
on the first Pareto front, namely the convex portion of the front. This provides 
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Figure 1: First 5 Pareto fronts on a toy example of 50 samples uniformly dis- 
tributed on [0, 1] X [0, 1]. The fronts create onion- like layers of the data. Samples 
on the first Pareto front indicated by a red '*' are minimizers under some linear 
combination. 



some intuition as to why the proposed method based on Pareto fronts is able to 
outperform methods using linear combinations. 

2.2 Related work 

Several machine learning methods utilizin g Pareto optimali t y hav e previously 



been proposed; an overview can be found in Jin and Sendhofi ( 20081) . These me- 



thods typically formulate machine learning problems as difficult multi-objective 
optimization problems. The goal is then to obtain a set of Pareto-optimal so- 
lutions using evolutionary algorithms or other methods for solving these multi- 
objective optimization problems. These methods differ from our use of Pareto 
optimality in that they use only the first Pareto front, while we use multiple 
Pareto fronts. Furthermore, we consider Pareto fronts created from a finite set 
of items, so we do not need to employ sophisticated methods in order to generate 

th ese fronts. 

Hero and Fleurv ( 20041) introduced a method for gene ranking and filtering 



using Pareto fronts that is closely related to our approach. The method es- 
tablishes a ranking of genes according to the estimated probability that each 
gene is Pareto-optimal. It differs from the method proposed in this paper in 
that the Pareto fronts are computed on the data samples themselves, i.e. the 
genes, rather than on dyads associated with the data samples. We introduce 
the concept of a dyad in Section 13.11 

Another related area is multi-view learning ( Blum and MitcheUl . 19981 : ISindhwani et al.l 



200i), which involves learning from data where observations are represented by 



multiple sets of features, commonly referred to as "views" . In such case, train- 
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ing in one view helps to improve learning in another view. The problem of 
view disagreement, where samples take different cla sses in different views, has 



recently been investigated (jChristoudias et al.l . 120081 ). The views are similar to 



criteria in our problem setting. However, in our setting, different criteria may 
be orthogonal and could even give contradictory information; hence there may 
be severe view disagreement. As a result, training in one view could actually 
worsen performance in another view, so the problem we consider differs from 
multi-view learning. 

Fin ally, many other anomaly dete ction methods h a ve pre viously been pro- 
posed. Hodge and Austin ( 2004 ) and IChandola et aL ( 2009[ ) both provide ex- 
tensive surveys of different anomaly detection methods and applications. Near- 
est neighbor-based metho ds are closely related to the proposed PDA approach. 
Bvers and Rafterv (Il998l) proposed to use the distance between a s ample and its 



fcth-ne arest neighbor as the ano maly score for the sample; similarly. I Angiulli and Pizzuti 
( 2002f ) and Eskin et al. ( 20021 ) proposed to the use the sum of the distances 
between a sample and its fc nearest neighbors. iBreunig et al. used an 

anomaly sco re based on the local density of the fc nearest neighbors of a sample. 
Herd ( 20061 ) introduced a non-parametric adaptive anomaly detection method 



using geometric ent ropy minimization, which is based on random fc-point mini- 
mal spanning trees. Zhao and Saligrama ( 20091 ) proposed an anomaly detection 
algorithm k-LPE using local p-value estimation (LPE) based on a fc-nearest 
neighbor graph. This algorithm is provably optimal for a specified level when 
the test samples are drawn from a mixture of the nominal density and a uniform 
density. 

All of the aforementioned methods are designed for single-criteria anomaly 
detection. In the multi-criteria case, the single-criteria algorithms must be exe- 
cuted multiple times with different weights, unlike the PDA anomaly detection 
algorithm that we propose. 



3 Pareto depth analysis 
3.1 Problem description 

Assume that a training set A'yv = {Xi, . . . ,Xn} of d-dimensional vectors Xi is 
available. Given a new vector X , the objective of anomaly detection is to declare 
X to be an anomaly if X is significantly different from samples in Ajv . Suppose 
that K > 1 different evaluation criteria are given. Denote these criteria by 
Ci , . . . , Ck ■ Each criterion is associated with a specific measure for computing 
dissimilarities. Denote the dissimilarity between Xi and Xj computed using the 
measure corresponding to criterion Ck by dk{i,j). 

We define a dyad by Aj = j), . . . ,dfe(i, j)]"^ e M+,i G {!,..., iV},j G 
{1, . . . , N} \ i. Each dyad Dij corresponds to a connection between vectors Xi 
and Xj. Therefore, there are in total (^) different dyads. For convenience, 
denote the set of all dyads by V and the space of all dyads M.^ by D. By 
the definition of strict dominance in Section [53] a dyad D strictly dominates 
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another dyad D* if D{i) < D*{i) for all i £ {l,...,K} and D{i) < D*{i) for 
some i, where D{i) denotes the ith component of D. Recall that each dyad 
corresponds to a connection between two samples. Therefore, if Dij >- Dik, the 
dissimilarity between Xi and Xj is no greater than that between Xi and X^ 
under any combination of the K criteria. The first Pareto front Fi corresponds 
to the set of dyads from T) that are not strictly dominated by any other dyads 
from T). The second Pareto front J-2 corresponds to the set of dyads from T)\IFi 
that are not strictly dominated by any other dyads from V\Ti, and so on, as 
defined in Section Recall that we refer to J^i as a deeper front than if 
i > j. 

3.2 Connection to anomaly detection 

If a dyad Dij is located in a shallow front, then the dissimilarity between the two 
samples Xi and Xj corresponding to Dij would be small under some combina- 
tion of the K criteria. For each sample X„, there are iV — 1 dyads corresponding 
to its connections with the other — 1 samples. Define the set of — 1 dyads 
associated with X„ by 2?". If most dyads in T>" are located at shallow Pareto 
fronts, then the dissimilarities between X„ and the other A^ — 1 samples are small 
under certain combinations of the criteria. Thus, X„ is likely to be a nominal 
sample. This is the basic idea of the proposed anomaly detection method using 
PDA and will play an important role in the following section. 

We construct Pareto fronts J^i, . . . ,J-m of the dyads from the training set, 
where the total number of fronts M is the required number of fronts such that 
each dyad is a member of a front. When a test sample X is obtained, we 
create new dyads corresponding to connections between X and training samples. 
Similar to many other anomaly detection methods, we connect each test sample 
to its k nearest neighbors. Since we are dealing with multiple criteria, we connect 
each test sample to its k nearest neighbors with respect to each criterion, and 
we do so using different values of k for each criterion, which we denote by 
fci,...,fcx. We delay the discussion of how to choose these fc^'s to Section 
14.21 We create s = J^iLi^^i ^^^w dyads, which we denote by the set 2?"°'^ = 
|£)ncw^ ^now^ _ _ ^ Corresponding to the connections between X and the 

union of the ki nearest neighbors in each criterion In other words, we create 
a dyad between X and Xj if Xj is among the ki nearest neighbors of X in any 
criterion i. We say that Df™ is below a front J^i if 

1?^^ >- Di for some Di e Fu 
i.e. strictly dominates at least a single dyad in J^;. Define the depth of 

e,: ^ min{/ 1 L>5^™ is below J",}. (1) 

^If a training sample is one of the ki nearest neighbors in m criteria, then m copies of the 
dyad corresponding to the connection between the test sample and the training sample are 
created. 
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Figure 2: Top: 40 training samples (blue circles) generated from bivariate mix- 
ture density and two test samples ('*' and '+')• Left: All dyads for the training 
samples (blue dots) along with first 20 Pareto fronts (green lines). Dyads asso- 
ciated with test point '*' (red circles) concentrate around shallow fronts. Right: 
Dyads associated with test point '+' (red circles) concentrate around deep fronts. 



Therefore if is large, then D"°'^ is around deep fronts, so the distance be- 
tween X and the corresponding training sample is large under all combinations 
of the K criteria. If is small, then Df°™ is around shallow fronts, so the 
distance between X and the corresponding training sample is small under some 
combination of the K criteria. 

To give a concrete example. Figure [5] shows the positions of 40 training 
samples and 2 test samples in along with the first 20 Pareto fronts using two 
criteria, the absolute value of the differences in the x and y dimensions. For 
demonstrative purposes, the test samples are connected to all of the training 
samples, i.e. ki — k2 = N. From Figure [U the is likely to be anomalous 
under any weighting of the two criteria while the '*' is likely to be nominal. 
Notice that the distributions of dyads in D associated with these two test samples 
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are quite different. The difference in these distributions is used in PDA to detect 
anomahes through the creation of an anomaly score. 

4 Multi-criteria anomaly detection 

Up to now, we have shown that the distribution of dyads in D with respect 
to the Pareto fronts can be very different for an anomalous sample compared 
to a nominal one. We formalize this idea to develop a multi-criteria anomaly 
detection algorithm. 

4.1 Anomaly score 

We first define an anomaly score for each test sample. We start with the ob- 
servation from Figure [2] that dyads corresponding to a nominal test sample are 
typically located near shallower fronts than dyads corresponding to an anoma- 
lous test sample. Therefore, the cumulative distribution function (cdf) of 
corresponding to an anomalous sample should be quite different from that of Cj 
corresponding to a nominal sample. Given 2?"'='"^ let z(i) denote the fraction of 
dyads in D"''^ below J", for i = 1, . . . , M. That is, 

ll^ncw g pncw | jjncw -g ^^^^^ jr.|| 



From ([T]), it can be seen that z is the cdf of the e^'s for a particular test sample. 
A plot of z as a function of front depth for the two test samples in Figure [2] is 
shown in Figure [31 Notice that the two cdf's do indeed differ significantly. The 
cdf for the anomalous sample grows more slowly as a function of the front 
depth because the dyads corresponding to the '+' are around deeper fronts as 
compared to the dyads corresponding to the 

For each test sample X, we define an anomaly score v{X) to capture the 
information contained in the cdf corresponding to X. We choose the anomaly 
score to be the area above the cdf of corresponding to X. That is, 

M 

v{X)^Y.^l~z{i)). 

i=l 

A benefit of this choice for anomaly score is that it can be calculated without 
first computing the cdf of by noting that it is simply the mean of the e^'s. This 
is easily shown using the identity E[ei] = Xi^i Pr(ei > j). Thus the anomaly 
score can be easily computed and compared to the decision threshold using the 
test 

v{X)^-Y,e,^o. 
s ^-^ 

i=l 

X is declared to be an anomaly if v{X) > a, where a is the decision threshold. 
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Figure 3: Cumulative distribution functions (cdf 's) of e^'s corresponding to test 
samples and '+' from Figure [2] The cdf for the grows more slowly because 
the dyads corresponding to the '+' are around deeper fronts. 



4.2 Selection of parameters 

Recall from Section [3^ that we connect each test sample X to a training sample 
Xj if Xj is one of the ki nearest neighbors of X in terms of the dissimilarity 
measure defined by criterion i for any i G {1,...,-?^}. We now discuss how 
these parameters fci, . . . , kx can be selected. For simplicity, first assume that 
there is only one criterion, so that a single parameter k is to be selected. PDA 
is able to detect an anomaly if the distribution of its dyads with respect to the 
Pareto fronts differs from that of a nominal sample. More formally, the area 
above the cdf of corresponding to an anomalous sample must be higher than 
that of a nominal sample. If k is chosen too small, this may not be the case, 
especially if there are training samples present near an anomalous sample, in 
which case, the dyads corresponding to the anomalous sample may reside near 
shallow fronts much like a nominal sample. On the other hand, if k is chosen 
too large, many dyads may correspond to connections to training samples that 
are far away, even if the test sample is nominal, which also makes the cdf 's of 
nominal and anomalous samples more similar. 

We propose to use the properties of /c-nearest neighbor graphs (/cNNGs) 
constructed on the training samples to select the number of training samples to 
connect to each test sample. We construct symmetric /cNNGs, i.e. we connect 
samples i and j if i is one of the k nearest neighbors of j or j is one of the k 
nearest neighbors of i. We begin with k — 1 and increase k until the fcNNG of the 
training samples is connected, i.e. there is only a single connected component. 
By forcing the fcNNG to be connected, we ensure that there are no isolated 
regions of training samples. Such isolated regions could possibly lead to dyads 
corresponding to anomalous samples residing near shallow fronts like nominal 
samples, which is undesirable. By keeping k small while retaining a connected 
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fcNNG, we are trying to avoid the problem of having too many dyads so that 
even a nominal sample may have many dyads located near deep fronts. This 
method of choosing k to retain connectivity has been used as a heuristic in othe r 
machine learning problems, such as spectral clustering (|von Luxburd . l2007l ). 
Note that by requiring the fcNNG to be connected, we are effectively imposing 
a one class assumption. If the nominal distribution consists of multiple well- 
separated classes, such an approach may not work well. 

Now let's return to the situation PDA was designed for, with K different 
criteria. For each criterion i, we construct a fc^NNGs using the correspond- 
ing dissimilarity measure and increase fci until the fc^NNG is connected. We 
then connect each test sample to s = X^i^i training samples. Note that we 
are choosing each ki independent of the other criteria, which is likely subop- 
timal. In principle, an approach that chooses the fc^'s jointly could perform 
better; however, such an approach would add to the complexity. Note that 
we choose separate ki's for each criterion, which we find is necessary to obtain 
good performance when different dissimilarities have varying scales and proper- 
ties. There are, however, pathological examples where the independent approach 
could choose fc^'s poorly, such as the well-known example of two moons. These 
examples typically involve multiple classes, which break the one-class assump- 
tion we effectively impose by requiring a connected fc^NNG. How to choose the 
fci's in a multi-class environment is beyond the scope of this paper and is an 
area for future work. 



4.3 Algorithm 

Pseudocode for training the PDA anomaly detector is shown in Algorithm [TJ It 
involves creating (^) dyads corresponding to all pairs of training samples then 
selecting the values of the parameters fci , . . . , fcx that determine how many 
training samples each test sample is connected to. Computing all pairwise dis- 
similarities in each criterion requires 0{dKN^) floating-point operations (flops), 
where d denotes the number of dimensions involved in computing a dissimilar- 
ity. Construction of the Parcto fro nts is done using n on-dominated sorting. We 
use the fast non-dominated sort of iDeb et al. ( 20001 ) that constructs all of the 
Pareto fronts using 0{KN'^) comparisons in the worst case where each front 
consists of a single dyad. 

Pseudocode for detecting whether a test sample is anomalous is shown in 
Algorithm [2] It involves creating dyads between the test sample and the ki 
nearest training samples in criterion I, which requires 0{dKN) flops. For each 
dyad Df™, we need to calculate the depth e^. This involves comparing the test 
dyad with training dyads on multiple fronts until we find a training dyad that is 
dominated by the test dyad, is the front that this training dyad is a part of. 
Using a binary search to select the front and another binary search to select the 
training dyads within the front to compare to, we need to make O(A'log^iV) 
comparisons in the worst case to compute e,;. The anomaly score is computed by 
taking the mean of the s e,; 's corresponding to the test sample; the score is then 
compared against a threshold a to determine whether the sample is anomalous. 
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Algorithm 1 PDA training phase. 
1: for I = 1 ^ K do 

2: Calculate pairwise dissimilarities di(i,j) between all training samples Xi 

and Xj 
3; fcz <- 
4: repeat 

5: ki^ki + 1 

6: Create fcjNNG using di(i,j)'s 
7: until fcjNNG is connected 
8: end for 

9: Create dyads Dij — [di{i,j), . . . ,(i_ff(j, j)] for all training samples 
10: Construct Pareto fronts on set of all dyads until each dyad is in a front 



Algorithm 2 PDA testing phase. 
1: nb ^ [] {empty list} 
2: for I = 1 ^ K do 

3: Calculate pairwise dissimilarities between test sample X and all training 

samples in criterion I 
4: nbi 4- ki nearest neighbors of X 
5: nb <— [nb, nbi] {append neighbors to list} 
6: end for 

7: Create s new dyads Df™ between X and training samples in nb 
8: for i = 1 — > s do 
9: Calculate depth e, of Df"" 
10: end for 

11: v{X)^{l/s)j:Ue, 

12: Declare X an anomaly if v{X) > a 



To handle multiple criteria, other anomaly detection methods, such as the 
ones mentioned in Section 12. 2[ need to be re-executed multiple times using 
different linear combinations of the K criteria. If a grid search is used for 
selection of the weights in the linear combination, then the required computation 
time would be exponential in K. Such an approach presents a computational 
problem unless K is very small. On the other hand, both the training time and 
the time required to test a new sample using PDA are linear in the number 
of criteria K. The PDA approach is clearly superior for large K; hence it is 
applicable to complex data sets where a large number of criteria must be used 
in order to achieve good performance. 

5 Experiments 

We compare the PDA method with several other nearest neighbor-based single- 
criterion anomaly detection algorithms including /cth-nearest neighbor distance 
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Table 1: Mean AUG (± standard error) comparison of different methods for the 
four-criteria simulation experiment. Best AUG (within one standard error) is 
in bold. 



Method 


25th 


AUG by percentile 
50th 75th 


100th 


PDA 




0.946 ± 0.004 




kNN 


0.804 ± 0.004 


0.856 ± 0.005 0.881 ± 0.005 


0.927 ± 0.005 


kNN sum 


0.809 ± 0.004 


0.859 ± 0.005 0.882 ± 0.006 


0.919 ± 0.005 


k-LPE 


0.802 ± 0.005 


0.854 ± 0.005 0.879 ± 0.005 


0.926 ± 0.005 


LOF 


0.791 ± 0.004 


0.845 ± 0.005 0.875 ± 0.005 


0.933 ± 0.004 



(kNN) (Ib^ ers and Raf terv'. ' 19 981), the sum of the distances to the k nearest 
neighbors (kNN sum) (Eskin et all . 120021: lAngiulh and Pizzuti l2002h . the k- 



Localized p- value Estimat or (k-LPE) ( Zhao and SaligramaL l2009t ) , and the Lo- 
cal Outlier Factor (LOF) (|Breunig et all . l2000l ). For these methods, we use lin- 
ear combinations of the criteria with different weights to compare performance 
with PDA. The linear combinations are selected using grid search. 



5.1 Simulated data with four criteria 

First we present an experiment on a simulated data set. The nominal distribu- 
tion is given by the uniform distribution on the hypercube [0,1]''. The anoma- 
lous samples are located just outside of this hypercube. There are four classes of 
anomalous distributions. Each class differs from the nominal distribution in one 
of the four dimensions; the distribution in the anomalous dimension is uniform 
on [1, 1.1]. We draw 300 training samples from the nominal distribution followed 
by 100 test samples from a mixture of the nominal and anomalous distributions 
with a 0.05 probability of selecting any particular anomalous distribution. The 
four criteria for this experiment correspond to the squared differences in each 
dimension. If the criteria are combined using linear combinations, the combined 
dissimilarity measure reduces to weighted squared Euclidean distanceH. 

The different methods are evaluated using the receiver operating character- 
istic (ROG) curve and the area under the curve (AUG). The mean AUCs (with 
standard errors) over 50 simulation runs are shown in Table [1] A grid of six 
points between and 1 in each criterion, corresponding to 6* = 1296 different 
sets of weights, is used to select linear combinations for all methods besides 
PDA. We summarize the results with quantile statistics. PDA is the best per- 
former; the other four methods all perform roughly equally. Notice that the 
performance of the single-criterion methods are indeed quite sensitive to the 
choice of weights. 

■^We obtain similar performance using the absolute differences in each dimension, for which 
linear combinations correspond to weighted £i distance. 
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5.2 Pedestrian trajectories 



We no w present an experiment on real data. This data set, provided bv lMaiecka 
( 2009| ). contains thousands of pedestrians' trajectories in an open area monitored 
by a video camera. E ach trajectory is approxi mated by a cubic spUne curve with 
seven control points ( Sillito and Fishen . 2008[) . We represent a trajectory by 



T ■ 



Xl X2 ... Xl 

yi y2 ■■■ yi 



where / denotes the number of time samples captured for the trajectory, and 
[xt, yt] denote a pedestrian's position at time step t. 

We use two criteria for computing the dissimilarity between trajectories. 
The first method is to compute the dissimilarity in walking speed. We com- 
pute the instantaneous speed at all time steps along each trajectory by fi- 
nite differencing, i.e. the speed of trajectory T at time step t is given by 
a/ {xt — xt^iy + {yt — Ut-iY- A histogram of speeds for each trajectory is ob- 
tained in this manner. We take the dissimilarity between two trajectories to 
be the squared Euclidean distance between their speed histograms. The sec- 
ond method is to compute the dissimilarity in shape. For each trajectory, we 
select 100 points, uniformly positioned along the trajectory. The dissimilarity 
between two trajectories T and T' is then given by the sum of squared Euclidean 
distances between the positions of T and T' over all 100 points. 

Because the scales and properties of these two dissimilarity measures are 
quite different, there is no simple way to normalize these two obtained dissim- 
ilarity matrices or to choose a "good" weight. Therefore if a single-criterion 
anomaly detection algorithm is executed separately for each dissimilarity mea- 
sure, it may miss some abnormal trajectories that are different from normal 
trajectories under a combination of both criteria but not very different under 
either of them individually. One would need to run such an algorithm with many 
linear combinations of the two criteria to detect a greater range of anomalies. 
However the proposed PDA method is able to do this without requiring selection 
of weights. 

The training sample for this experiment consists of 500 trajectories, and 
the test sample consists of 200 trajectories. Figure |4] shows some abnormal 
trajectories and nominal trajectories detected using PDA. Recall that anoma- 
lous trajectories could have anomalous speeds or shapes (or both), so some 
anomalous trajectories in Figure |4] may not look anomalous by shape alone. We 
select parameters in PDA using the method described in Section W?^ resulting 
in [fci = 3, fc2 — 6] for the speed and shape criteria, respectively. The AUCs 
of the PDA method using different [fci,fc2] pairs are shown in Figure [S] The 
selected [fci, fc2] pair is very close to the optimal [fci, ^2] pair [3, 7], and the AUG 
is also very close. Table [2] shows the performance of PDA as compared to the 
other algorithms using 100 different choices of weights. Notice that PDA has 
larger AUG than the other methods under all choices of weights for the two 
criteria. For a more detailed comparison, the ROG curve for PDA and the at- 
tainable region for k-LPE (the region between the ROC curves corresponding 
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Figure 4: Left: Some abnormal trajectories detected by PDA method. Right: 
Trajectories with relatively low anomaly scores. 
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Figure 5: AUG for different choices of [A:i,fc2]. The automatically selected pa- 
rameters [ki ^ 3, krt ~ 6] are very close to the optimal parameters [ki — 3, k2 = 
7]. 



to the best and worst AUCs) is shown in Figure [6] Notice that k-LPE performs 
slightly better at low false positive rate when the best weights are used, but 
PDA performs better in all other situations, resulting in higher AUG. 



6 Conclusion 

In this paper we proposed a new multi-criteria anomaly detection method. The 
proposed method uses Pareto depth analysis to compute the anomaly score of 
a test sample by examining the Pareto front depths of dyads corresponding 
to the test sample. Dyads corresponding to an anomalous sample tended to 
be located at deeper fronts compared to dyads corresponding to a nominal 
sample. Instead of choosing a specific weighting or performing a grid search 
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Table 2: AUG comparison of different methods for the trajectory experiment. 
Best AUG is in bold. 



Method AUG by percentile 

^ ° 25th 50th 75th 100th 



PDA 




0.915 




kNN 


0.869 


0.883 


0.900 


0.906 


kNN Sum 


0.878 


0.894 


0.908 


0.911 


k-LPE 


0.874 


0.893 


0.904 


0.908 


LOF 


0.772 


0.800 


0.810 


0.833 




False positive rate 



Figure 6: ROG curves for PDA and attainable region for k-LPE over 100 choices 
of weights. PDA outperforms k-LPE even under the best choice of weights. 



on the weights for different dissimilarity measures, the proposed method can 
efficiently detect anomalies under multiple combinations at once. The proposed 
method scales linearly in the number of criteria, which makes it applicable to 
complicated data sets where a large number of criteria may be necessary to 
identify anomalies. 

References 

F. AngiuUi and G. Pizzuti (2002). Fast outlier detection in high dimensional 
spaces. In Proceedings of the 6th European Conference on Principles of Data 
Mining and Knowledge Discovery. 

A. Blum and T. Mitchell (1998). Gombining labeled and unlabeled data with 



15 



co-training. In Proceedings of the 11th Annual Conference on Computational 

Learning Theory. 

M. M. Breunig, H.-P. Kriegel, R. T. Ng, and J. Sander (2000). LOF: Identifying 
density-based local outliers. In Proceedings of the ACM SIGMOD International 
Conference on Management of Data. 

S. Byers and A. E. Raftery (1998). Nearest-neighbor chittcr removal for esti- 
mating features in spatial point processes. Journal of the American Statistical 

Association 93{U2):r:>77 584. 

V. Chandola, A. Banerjee, and V. Kumar (2009). Anomaly detection: A survey. 
ACM Computing Surveys 41(3):l-58. 

C. M. Christoudias, R. Urtasun, and T. Darrell (2008). Multi-view learning in 

the presence of view disagreement. In Proceedings of the Conference on Uncer- 
tainty in Artificial Intelligence. 

K. Deb, S. Agrawal, A. Pratap, and T. Meyarivan (2000). A fast elitist non- 
dominated sorting genetic algorithm for multi-objective optimization: NSGA-II. 
In Proceedings of the 6th International Conference on Parallel Problem Solving 
from Nature. 

E. Eskin, A. Arnold, M. Prcrau, L. Portnoy, and S. Stolfo (2002). A geometric 
framework for unsupervised anomaly detection: Detecting intrusions in unla- 
beled data. In Applications of Data Mining in Computer Security. Norwell, MA: 
Kluwer. 

A. O. Hero III and G. Fleury (2004). Parcto-optimal methods for gene ranking. 
The Journal of VLSI Signal Processing 38(3):259-275. 

A. O. Hero HI (2006). Geometric entropy minimization (GEM) for anomaly de- 
tection and localization. In Advances in Neural Information Processing Systems 
19. 

V. J. Hodge and J. Austin (2004). A survey of outlier detection methodologies. 
Artificial Intelligence Review 22(2):85-126. 

Y. Jin and B. Sendhoff (2008). Pareto-bascd multiobjectivc machine learning: 
An overview and case studies. IEEE Transactions on Systems, Man, and Cy- 
bernetics, Part C: Applications and Reviews 38(3):397-415. 

B. Majecka (2009). Statistical models of pedestrian behaviour in the Forum. 
Master's thesis, University of Edinburgh. 

R. R. Sillito and R. B. Fisher (2008). Semi-supervised learning for anomalous 
trajectory detection. In Proceedings of the 19th British Machine Vision Confer- 
ence. 



16 



V. Sindhwani, P. Niyogi, and M. Belkin (2005). A co-rcgularization approach to 
semi-supervised learning with multiple views. In Proceedings of the Workshop 
on Learning with Multiple Views, 22nd International Conference on Machine 
Learning. 

U. von Luxburg (2007). A tutorial on spectral clustering. Statistics and Com- 
puting 17 (4) -.395 416. 

M. Zhao and V. Saligrama (2009). Anomaly detection with score functions based 
on nearest neighbor graphs. In Advances in Neural Information Processing Sys- 
tems 22. 



17 



